Exercise 1

a)

# We will work with companies like Apple, Microsoft and Intel
symbols <- c("AAPL", "MSFT", "INTC")
# Using the functions from functions_Assignment_3.R to get stock data
df_all <- get_finance_data(symbols, "2000-01-01", Sys.Date())

head(df_all)
## # A tibble: 6 × 5
##   date       symbol  open   low adjusted
##   <date>     <fct>  <dbl> <dbl>    <dbl>
## 1 2000-01-03 AAPL   0.936 0.908    0.840
## 2 2000-01-04 AAPL   0.967 0.903    0.769
## 3 2000-01-05 AAPL   0.926 0.920    0.781
## 4 2000-01-06 AAPL   0.948 0.848    0.713
## 5 2000-01-07 AAPL   0.862 0.853    0.747
## 6 2000-01-10 AAPL   0.911 0.846    0.734
monthly_data <- df_all %>%
  # Creates new column so every row knows which month it belongs to
  mutate(year_month = floor_date(date, "month")) %>%
  # Group by stock and month 
  group_by(symbol, year_month) %>%
  # Selects last day of each month
  slice_max(date, n = 1, with_ties = FALSE) %>%
  # Keeps only the columns date, symbol, and adjusted
  dplyr::select(date, symbol, adjusted) %>%
  ungroup() %>%
  # Sort by stock alphabetically and date 
  arrange(symbol, date)

head(monthly_data)
## # A tibble: 6 × 4
##   year_month date       symbol adjusted
##   <date>     <date>     <fct>     <dbl>
## 1 2000-01-01 2000-01-31 AAPL      0.779
## 2 2000-02-01 2000-02-29 AAPL      0.860
## 3 2000-03-01 2000-03-31 AAPL      1.02 
## 4 2000-04-01 2000-04-28 AAPL      0.931
## 5 2000-05-01 2000-05-31 AAPL      0.630
## 6 2000-06-01 2000-06-30 AAPL      0.786
monthly_prices <- monthly_data %>%
  # Select date, stock and adjusted price, removing year_month column
  dplyr::select(date, symbol, adjusted) %>%
  # Transforms the data from long format to wide format:
  pivot_wider(names_from = symbol, values_from = adjusted) %>%
  # Sorts rows
  arrange(date)

head(monthly_prices)
## # A tibble: 6 × 4
##   date        AAPL  INTC  MSFT
##   <date>     <dbl> <dbl> <dbl>
## 1 2000-01-31 0.779  28.1  29.9
## 2 2000-02-29 0.860  32.1  27.3
## 3 2000-03-31 1.02   37.5  32.5
## 4 2000-04-28 0.931  36.0  21.3
## 5 2000-05-31 0.630  35.4  19.1
## 6 2000-06-30 0.786  38.0  24.5
summary(monthly_prices[, -1])
##       AAPL               INTC             MSFT       
##  Min.   :  0.2122   Min.   : 7.941   Min.   : 11.90  
##  1st Qu.:  2.1356   1st Qu.:14.160   1st Qu.: 19.16  
##  Median : 16.1177   Median :19.377   Median : 26.08  
##  Mean   : 48.2217   Mean   :24.407   Mean   : 98.76  
##  3rd Qu.: 51.4990   3rd Qu.:31.996   3rd Qu.:124.76  
##  Max.   :270.3700   Max.   :58.209   Max.   :532.62
# Fill missing values using last observation carried forward (LOCF) due to non-trading days
monthly_prices <- monthly_prices %>%
  mutate(across(everything(), ~na.locf(.x, na.rm = FALSE)))

# Confirm that NAs are gone
colSums(is.na(monthly_prices))
## date AAPL INTC MSFT 
##    0    0    0    0

Comments on summarized stock data:

Microsoft has the highest mean and median stock prices, indicating a consistently strong market position. Intel shows the lowest mean, while Apple has the lowest median, suggesting that Apple’s stock remained relatively low for much of the early 2000s before its rapid growth phase. The large gap between mean and median for both Microsoft and Apple reflects significant price increases in recent years, consistent with strong long-term growth in the tech sector.

b)

historical_events <- data.frame(
  date = as.Date(c("2001-09-11", "2008-10-01", "2020-03-15",
                   "2020-03-23", "2022-02-24", "2007-01-09",
                   "2015-07-29", "2006-11-02", "2008-09-23",
                   "2010-04-03", "2001-10-25", "2009-10-22" )),
  event = c("9/11 Attacks", "Financial Crisis", "COVID-19 Pandemic", 
            "COVID Market Bottom", "Russia-Ukraine War", "iPhone Announced",
            "Windows 10 Launch", "Intel Core 2 Duo Launch", "Intel Core i7 Launch",
            "iPad Launch", "Windows XP release", "Windows 7 launch"),
  y_pos = c(400, 15, 350, 20, 450, 25, 300, 30, 12, 380, 320, 18),
  period = c("2000s", "2008", "2020", "2020", "2020s", "2007", 
             "2010s", "2006", "2008", "2010", "2001", "2009")
)
monthly_data %>%
  ggplot(aes(x = date, y = adjusted, color = symbol)) +
  geom_line(size = 1.2) +
  geom_vline(data = historical_events, aes(xintercept = date),
             color = "gray40", linetype = "dashed", alpha = 0.6) +
  geom_text_repel(data = historical_events, 
                  aes(x = date, y = 100, label = event),
                  angle = 90, hjust = 0, vjust = 0.5, size = 3.5, 
                  color = "black", fontface = "bold",
                  bg.color = "white", bg.r = 0.1,
                  box.padding = 0.5, point.padding = 0.5,
                  force = 2, max.overlaps = Inf) +
  scale_y_log10(labels = scales::dollar_format()) +
  scale_color_manual(values = c("AAPL" = "#1f77b4", "MSFT" = "#2ca02c", "INTC" = "#8c564b")) +
  labs(title = "Stock Price Evolution (Log Scale)", x = "Year", y = "Price (USD, log scale)") +
  theme_minimal() +
  theme(
    plot.title = element_text(size = 16, hjust = 0.5),
    axis.title = element_text(size = 14),
    axis.text = element_text(size = 12),
    legend.title = element_text(size = 14),
    legend.text = element_text(size = 12)
  )

9/11 Attacks (2001): The attacks caused short-term uncertainty in the stock market, and all three stocks dipped briefly. However, the impact was relatively modest and short-lived compared to later crises.

Intel product launches (Core 2 Duo, Core i7): Despite new product releases, Intel’s stock performance remained quite volatile. The Core i7 announcement in particular did not boost the stock, which actually declined afterward. It could be possible that they didnt meet the expecations.

Apple’s rise (2003 onward): From around 2003, Apple’s stock began a strong, sustained growth trend. The launches of the iPhone (2007) and iPad (2010) appear to have accelerated this upward momentum.

Financial Crisis (2008): All three stocks dropped sharply during the global financial crisis.

Microsoft’s rebound (Windows 10 launch): Following the release of Windows 10 (2015), Microsoft’s stock rose substantially.

COVID-19 pandemic (2020): Initially, stock prices fell sharply as markets reacted to uncertainty. However, the rebound was rapid — particularly for Apple and Microsoft. The demand for technology and digital servieces had still a strong demand during lockdowns.

Russia–Ukraine War (2022): The invasion triggered renewed volatility. Intel’s stock was hit hardest, possibly due to supply chain disruptions and its exposure to manufacturing and global chip markets, while Microsoft and Apple were less affected in the long run.

Over the past two decades, Apple has shown the strongest long-term growth, Microsoft recovered strongly after its mid-2010s while Intel has lagged behind. ### c)

log_returns <- monthly_data %>%
  arrange(symbol, date) %>%
  group_by(symbol) %>%
  mutate(log_return = log(adjusted / lag(adjusted))) %>%
  filter(!is.na(log_return) & !is.infinite(log_return)) %>%
  ungroup()

summary(log_returns)
##    year_month              date             symbol       adjusted       
##  Min.   :2000-02-01   Min.   :2000-02-29   AAPL:310   Min.   :  0.2122  
##  1st Qu.:2006-07-01   1st Qu.:2006-07-31   INTC:310   1st Qu.: 14.0898  
##  Median :2012-12-16   Median :2013-01-15   MSFT:310   Median : 21.2813  
##  Mean   :2012-12-15   Mean   :2013-01-13              Mean   : 57.2511  
##  3rd Qu.:2019-06-01   3rd Qu.:2019-06-28              3rd Qu.: 45.1086  
##  Max.   :2025-11-01   Max.   :2025-11-06              Max.   :532.6244  
##    log_return       
##  Min.   :-0.861415  
##  1st Qu.:-0.041497  
##  Median : 0.017140  
##  Mean   : 0.009611  
##  3rd Qu.: 0.065965  
##  Max.   : 0.374169
# Fill missing values using last observation carried forward (LOCF) due to non-trading days
log_returns <- log_returns %>%
  mutate(across(everything(), ~na.locf(.x, na.rm = FALSE)))

# Confirm that NAs are gone
colSums(is.na(log_returns))
## year_month       date     symbol   adjusted log_return 
##          0          0          0          0          0
# Time series plot of log returns for all companies
log_returns %>%
  ggplot(aes(x = date, y = log_return, color = symbol)) +
  geom_line(size = 0.6, alpha = 0.8) +
  geom_hline(yintercept = 0, linetype = "dashed", color = "gray50") +
  facet_wrap(~symbol, ncol = 1, scales = "free_y") +
  scale_color_manual(values = c("AAPL" = "#1f77b4", "MSFT" = "#2ca02c", "INTC" = "#8c564b")) +
  labs(title = "Monthly Log Returns Over Time", 
       x = "Year", 
       y = "Log Return",
       subtitle = "Separate plots for each company") +
  theme_minimal() +
  theme(
    plot.title = element_text(size = 16, hjust = 0.5),
    plot.subtitle = element_text(size = 12, hjust = 0.5),
    axis.title = element_text(size = 14),
    axis.text = element_text(size = 12),
    strip.text = element_text(size = 14, face = "bold"),
    legend.position = "none"
  )

The volatility looks very similar to all the companies by looking at the plot.

volatility_stats <- log_returns %>%
  group_by(symbol) %>%
  summarise(
    observations = n(),
    mean_return = mean(log_return, na.rm = TRUE),
    volatility = sd(log_return, na.rm = TRUE),
    annualized_vol = sd(log_return, na.rm = TRUE) * sqrt(12),
    min_return = min(log_return, na.rm = TRUE),
    max_return = max(log_return, na.rm = TRUE),
    .groups = 'drop'
  )

print(volatility_stats)
## # A tibble: 3 × 7
##   symbol observations mean_return volatility annualized_vol min_return
##   <fct>         <int>       <dbl>      <dbl>          <dbl>      <dbl>
## 1 AAPL            310    0.0189       0.115           0.398     -0.861
## 2 INTC            310    0.000908     0.106           0.369     -0.589
## 3 MSFT            310    0.00906      0.0789          0.273     -0.421
## # ℹ 1 more variable: max_return <dbl>

Although the three companies show visually similar volatility patterns over time, Apple has the highest measured volatility, followed by Intel and then Microsoft. This means Apple’s stock price fluctuates more strongly month to month, reflecting higher risk and sensitivity to market news, but also the potential for higher returns. Microsoft appears to be the most stable of the three, with smaller swings in its monthly log returns.

d)

# Density plot comparison
p <- log_returns %>%
  ggplot(aes(x = log_return, color = symbol)) +
  geom_density(alpha = 0.3, size = 1.2) +
  scale_color_manual(values = c("AAPL" = "#1f77b4", "MSFT" = "#2ca02c", "INTC" = "#8c564b")) +
  labs(title = "Density Comparison of Log Returns", x = "Log Return", y = "Density") +
  theme_minimal() +
  theme(
    plot.title = element_text(size = 16, hjust = 0.5),
    axis.title = element_text(size = 14),
    axis.text = element_text(size = 12),
    legend.title = element_text(size = 14),
    legend.text = element_text(size = 12)
  )

print(p)

The distributions of monthly log returns for all three companies are roughly centered around zero, which is expected since stock returns fluctuate around an average long-term growth rate. None of them are perfectly symmetrical — they show slight negative skewness, meaning large negative returns occur a bit more frequently than large positive ones.

Apple’s curve is flatter and wider, reflecting its higher volatility (larger spread of returns). Microsoft’s distribution is narrower and more peaked, indicating more stable returns and lower volatility. Intel lies between the two, with moderate spread.

These shapes align with the earlier time-series plots from part 1(c): Apple’s returns fluctuate more strongly over time, Microsoft’s are steadier, and Intel shows intermediate variability. The heavier tails seen here also correspond to the occasional sharp spikes (extreme returns) visible in the time-series plots. ### e)

log_returns_wide <- log_returns %>%
  dplyr::select(date, symbol, log_return) %>%
  pivot_wider(names_from = symbol, values_from = log_return) %>%
  arrange(date) %>%
  filter(!is.na(AAPL) & !is.na(MSFT) & !is.na(INTC))

cor_matrix <- log_returns_wide %>%
  dplyr::select(AAPL, MSFT, INTC) %>%
  cor(use = "complete.obs")

print(round(cor_matrix, 4))
##        AAPL   MSFT   INTC
## AAPL 1.0000 0.4508 0.5461
## MSFT 0.4508 1.0000 0.4524
## INTC 0.5461 0.4524 1.0000

The correlation matrix shows that all three stocks are positively correlated, meaning they tend to move in the same direction over time.

The strongest relationship is between Intel (INTC) and Apple (AAPL) with a correlation of 0.55, suggesting that their returns often rise or fall together, possibly due to shared exposure to the broader tech industry and global market trends.

Microsoft (MSFT) has a more moderate correlation with both Apple (0.45) and Intel (0.45), indicating some common movement but also a degree of independence in how its stock reacts to market conditions.

Overall, the correlations are moderate rather than high, meaning while the stocks are influenced by similar market factors, each company’s specific performance and product strategies also play an important role in their return patterns.

Exercise 2

We wanted to further explore Apple’s stock, which has the highest volatility and highest acceleration in price.

a)

apple_data <- monthly_data %>%
# Keep only Apple Stock
  filter(symbol == "AAPL") %>%
  arrange(date) %>%
  dplyr::select(date, adjusted)

head(apple_data)
## # A tibble: 6 × 2
##   date       adjusted
##   <date>        <dbl>
## 1 2000-01-31    0.779
## 2 2000-02-29    0.860
## 3 2000-03-31    1.02 
## 4 2000-04-28    0.931
## 5 2000-05-31    0.630
## 6 2000-06-30    0.786
apple_log_returns <- log_returns %>%
  filter(symbol == "AAPL") %>%
  arrange(date) %>%
  dplyr::select(date, log_return)

# Replace NA values by carrying forward the last non-NA value
apple_data$adjusted <- zoo::na.locf(apple_data$adjusted)
apple_log_returns$log_return <- zoo::na.locf(apple_log_returns$log_return)

# Convert to time series
apple_ts <- ts(apple_data$adjusted, start = c(2000, 2), frequency = 12)
apple_returns_ts <- ts(apple_log_returns$log_return, start = c(2000, 3), frequency = 12)

A stationary time series is one whose statistical properties do not change over time. In other words, the process generating the data is stable and behaves similarly in the past, present, and future. For Apple’s stock price to be stationary, its mean and variance should remain constant over time. However, from earlier observations, the stock price shows a clear upward trend, meaning the mean changes over time. Therefore, Apple’s stock price is not stationary.

# KPSS test for prices
kpss.test(apple_ts)
## 
##  KPSS Test for Level Stationarity
## 
## data:  apple_ts
## KPSS Level = 3.8147, Truncation lag parameter = 5, p-value = 0.01
nsdiffs(apple_ts)
## [1] 0
# KPSS test for log returns
kpss.test(apple_returns_ts)
## 
##  KPSS Test for Level Stationarity
## 
## data:  apple_returns_ts
## KPSS Level = 0.098696, Truncation lag parameter = 5, p-value = 0.1
nsdiffs(apple_returns_ts)
## [1] 0

The KPSS test for Apple’s stock prices gives a very low p-value (< 0.05), indicating that the price series is non-stationary — it has a trend and changing mean over time. In contrast, the KPSS test for log returns shows a high p-value (> 0.05), suggesting that the return series is stationary. This makes sense, since stock prices typically follow a random walk (non-stationary ), while returns fluctuate around a stable mean with relatively constant variance. Therefore, returns are more suitable for time series modeling and forecasting than raw prices.

b)

# STL decomposition for Apple prices
stl_price <- stl(apple_ts, s.window = "periodic")
plot(stl_price, main = "STL Decomposition: Apple Stock Prices")

# STL decomposition for Apple log returns
stl_returns <- stl(apple_returns_ts, s.window = "periodic")
plot(stl_returns, main = "STL Decomposition: Apple Log Returns")

The STL decomposition of Apple’s stock price shows an upward trend but no meaningful seasonal pattern. Since stock prices do not follow consistent seasonal cycles, the seasonal component detected by STL likely reflects random fluctuations rather than true seasonality. This limits the interpretability of STL for stock data, though the method still effectively highlights the long-term trend and residual volatility.

c)

# Autocorrelation and Partial Autocorrelation Analysis

# Short-term analysis (1 year = 12 lags)
par(mfrow = c(2, 2))

# ACF and PACF for Apple prices - short term
acf(apple_ts, lag.max = 12, main = "ACF: Apple Prices (1 year)")
pacf(apple_ts, lag.max = 12, main = "PACF: Apple Prices (1 year)")

# ACF and PACF for Apple log returns - short term
acf(apple_returns_ts, lag.max = 12, main = "ACF: Apple Log Returns (1 year)")
pacf(apple_returns_ts, lag.max = 12, main = "PACF: Apple Log Returns (1 year)")

ACF & PACF for Apple Prices (1 year analysis):

ACF: The autocorrelation values are very high and decay slowly, nearly 1 at lag 1 and still large after 12 lags. This means each month’s price is highly correlated with previous months’ prices — typical for a non-stationary series with a strong upward trend.

PACF: The PACF for Apple’s prices shows one strong spike at lag 1, indicating a high direct correlation with the previous month’s price. All higher-order lags fall within the confidence bounds, meaning that once the effect of the previous month is accounted for, earlier months have little additional influence.

ACF & PACF for Apple Log Returns (1 year analysis):

ACF: Only the first lag is near 1, and all others are small and within the blue confidence bounds. This means there’s little to no serial correlation, consistent with a stationary, random-like process.

PACF: Also shows only minor spikes, mostly within the bounds. No clear cutoff or strong structure appears.

Apple’s log returns are approximately stationary — their mean and variance are stable, and there’s no persistent correlation over time.

# Long-term analysis (20 years = 240 lags)
par(mfrow = c(2, 2))

# ACF and PACF for Apple prices - long term
acf(apple_ts, lag.max = 240, main = "ACF: Apple Prices (20 years)")
pacf(apple_ts, lag.max = 240, main = "PACF: Apple Prices (20 years)")

# ACF and PACF for Apple log returns - long term
acf(apple_returns_ts, lag.max = 240, main = "ACF: Apple Log Returns (20 years)")
pacf(apple_returns_ts, lag.max = 240, main = "PACF: Apple Log Returns (20 years)")

par(mfrow = c(1, 1))

ACF & PACF for Apple Prices (20 year analysis):

ACF: The ACF starts at 1.0 at lag 1 and decays very slowly over many lags, remaining positive for a long time. This indicates strong persistence — each month’s price is heavily correlated with prices from many months or even years earlier.

PACF: The PACF shows a single strong spike at lag 1, followed by very small and insignificant correlations beyond that point. This suggests that most of the price’s dependence structure can be explained by the previous month’s value — once that is accounted for, earlier lags add little explanatory power.

ACF & PACF for Apple Log Returns (1 year analysis):

ACF: The autocorrelation drops immediately to near zero after lag 1, and all subsequent lags are within the blue confidence bounds. This means there is no significant serial correlation in the returns — past monthly returns have almost no predictive power for future returns.

PACF: The PACF also shows no clear significant spikes, with nearly all values within the blue confidence bounds. This indicates that once any minor first-lag effect is removed, there’s no direct relationship between returns and earlier periods.

d)

The Granger causality test evaluates whether past values of one variable X contain useful information for predicting another variable Y. The null hypothesis states that X does not Granger-cause Y, meaning that the coefficients on lagged X terms in a VAR model are all zero. The test compares an autoregressive model of Y with a VAR model including both X and Y. A significant F-statistic indicates that adding lagged values of X improves the prediction of Y, so X Granger-causes Y. ### e)

# Granger causality tests for prices
price_data <- monthly_prices %>% dplyr::select(AAPL, MSFT, INTC) # Removes date column
var_model_prices <- VAR(price_data, p = 2) # Fits var model, p=2 is apparently common in finance
granger_test(var_model_prices) #Tests granger casuality relationship
## 
## Granger Causality Test (Multivariate)
## 
## F test based on VAR(2) model:
## ----------------------------------------
##                      F df1 df2     p    
## ----------------------------------------
##  ------------                           
##  AAPL <= MSFT    25.81   2 302 <.001 ***
##  AAPL <= INTC     5.28   2 302  .006 ** 
##  AAPL <= ALL     16.83   4 302 <.001 ***
##  ------------                           
##  MSFT <= AAPL     3.69   2 302  .026 *  
##  MSFT <= INTC     2.51   2 302  .083 .  
##  MSFT <= ALL      3.93   4 302  .004 ** 
##  ------------                           
##  INTC <= AAPL     0.60   2 302  .549    
##  INTC <= MSFT     0.11   2 302  .896    
##  INTC <= ALL      0.69   4 302  .601    
## ----------------------------------------

The results shows that Microsoft and Intel both Granger-cause Apple’s stock prices (p < 0.01), indicating that their past prices contain predictive information about Apple. Microsoft’s prices are weakly Granger-caused by Apple (p ≈ 0.045), while Intel’s prices show no significant causal links with either Apple or Microsoft. Overall, Apple appears most sensitive to the price dynamics of the other two firms, whereas Intel’s price movements are largely independent.

# Granger causality tests for log returns
returns_data <- log_returns_wide %>% dplyr::select(AAPL, MSFT, INTC)
var_model_returns <- VAR(returns_data, p = 2)
granger_test(var_model_returns)
## 
## Granger Causality Test (Multivariate)
## 
## F test based on VAR(2) model:
## ---------------------------------------
##                     F df1 df2     p    
## ---------------------------------------
##  ------------                          
##  AAPL <= MSFT    0.69   2 301  .503    
##  AAPL <= INTC    2.00   2 301  .137    
##  AAPL <= ALL     1.20   4 301  .311    
##  ------------                          
##  MSFT <= AAPL    2.05   2 301  .130    
##  MSFT <= INTC    0.44   2 301  .642    
##  MSFT <= ALL     1.04   4 301  .388    
##  ------------                          
##  INTC <= AAPL    1.30   2 301  .274    
##  INTC <= MSFT    0.07   2 301  .928    
##  INTC <= ALL     0.68   4 301  .609    
## ---------------------------------------

The Granger causality tests on log returns show no significant causal relationships between Apple, Microsoft, and Intel. All p-values exceed 0.05, indicating that past returns of one company do not help predict future returns of another.

Comparing the results:

When using price levels, several significant causal relationships appear. It could be that the result reflects the shared upward trends and long-term co-movements of the companies’ stock prices rather than genuine predictive power. When the same test is applied to log returns, all relationships become statistically insignificant, indicating that past returns do not help predict future returns across companies.

f)

Granger causality does not imply true causality. It only indicates that one variable’s past values improve the prediction of another variable’s current value. The relationship is predictive rather than structural — other hidden factors or common trends may drive both variables. Therefore, while Granger tests reveal useful information about temporal dependencies, they cannot establish real cause-and-effect relationships.

Exercise 3

a)

train-test split:

h <- 48
stopifnot(frequency(apple_returns_ts) == 12)

n <- length(apple_returns_ts)
if (n <= h) stop("For få observasjoner til å ha 4 års testsett.")

ts_time <- time(apple_returns_ts)

apple_train <- window(apple_returns_ts, end   = ts_time[n - h])
apple_test  <- window(apple_returns_ts, start = ts_time[n - h + 1])

b)

Give a brief explanation of the 4 baseline methods:

• average method - This method forecasts future values by calculating the mean of all past observations. • drift method - This method forecasts future values by extending a line from the first to the last observation, effectively capturing the overall trend. • naive method - This method forecasts future values by using the last observed value as the forecast for all future values. • seasonal naive method - This method forecasts future values by using the last observed value from the same season (e.g., month) as the forecast for all future values in that season.

What is the purpose of starting with these methods? What are the requirements for the residuals of time series models? - The purpose of starting with these baseline methods is to establish a simple benchmark for forecasting performance. - These methods are easy to implement and understand, providing a reference point against which more complex models can be compared. - The requirements for the residuals of time series models are that they should be uncorrelated (white noise), have a mean of zero, and exhibit constant variance (homoscedasticity).

c)

Train the 4 baseline models:

avg_fc   <- forecast::meanf(apple_train, h = h)
drift_fc <- forecast::rwf(apple_train, h = h, drift = TRUE)
naive_fc <- forecast::rwf(apple_train, h = h, drift = FALSE)
seasonal_naive_fc <- forecast::snaive(apple_train, h = h)

Plot forecast

autoplot(apple_train, series = "Training data") +
  autolayer(apple_test, series = "Test data") +
  autolayer(avg_fc$mean,   series = "Average method") +
  autolayer(drift_fc$mean, series = "Drift method") +
  autolayer(naive_fc$mean, series = "Naive method") +
  autolayer(seasonal_naive_fc$mean, series = "Seasonal naive method") +
  ggtitle("Apple Monthly Log Returns – Baseline Forecasts") +
  xlab("Year") + ylab("Log Return") +
  theme_minimal() +
  scale_color_manual(values = c(
    "Training data" = "black",
    "Test data" = "gray50",
    "Average method" = "#1f77b4",
    "Drift method" = "#ff7f0e",
    "Naive method" = "#2ca02c",
    "Seasonal naive method" = "#d62728"
  )) +
  theme(
    plot.title = element_text(size = 15, hjust = 0.5),
    axis.title = element_text(size = 13),
    legend.title = element_blank(),
    legend.position = "bottom"
  )

The base models gives expected results. With the average method, the forecast is a flat line at the mean of the training data. With the drift method, the forecast extends the trend observed in the training data. With the naive method, the forecast is a flat line at the last observed value of the training data. With the seasonal naive method, the forecast repeats the last observed seasonal pattern, resulting in large fluctuations that closely follow recent values.

compute RMSE and MAE

avg_rmse <- sqrt(mean((apple_test - avg_fc$mean)^2))
avg_mae <- mean(abs(apple_test - avg_fc$mean))

drift_rmse <- sqrt(mean((apple_test - drift_fc$mean)^2))
drift_mae <- mean(abs(apple_test - drift_fc$mean))

naive_rmse <- sqrt(mean((apple_test - naive_fc$mean)^2))
naive_mae <- mean(abs(apple_test - naive_fc$mean))

seasonal_naive_rmse <- sqrt(mean((apple_test - seasonal_naive_fc$mean)^2))
seasonal_naive_mae <- mean(abs(apple_test - seasonal_naive_fc$mean))

RMSE and MAE

results <- data.frame(
  method = c("Average", "Drift", "Naive", "Seasonal Naive"),
  RMSE   = round(c(avg_rmse, drift_rmse, naive_rmse, seasonal_naive_rmse), 6),
  MAE    = round(c(avg_mae,  drift_mae,  naive_mae,  seasonal_naive_mae),  6),
  stringsAsFactors = FALSE
)

header <- sprintf("%-15s %12s %12s", "method", "RMSE", "MAE")
rows <- apply(results, 1, function(row)
  sprintf("%-15s %12.6f %12.6f", row[1], as.numeric(row[2]), as.numeric(row[3]))
)
cat(header, "\n", paste(rows, collapse = "\n"), "\n", sep = "")
## method                  RMSE          MAE
## Average             0.071625     0.059019
## Drift               0.114078     0.095047
## Naive               0.114061     0.095030
## Seasonal Naive      0.089963     0.072704

method RMSE MAE

Average 0.071619 0.058998 Drift 0.114059 0.095026 Naive 0.114042 0.095009 Seasonal Naive 0.089940 0.072683

Comment: Among the baseline models, the Average method performs best with the lowest RMSE and MAE, while the Drift and Naive methods perform similarly and considerably worse. The Seasonal Naive method shows slightly better accuracy than Drift and Naive but still underperforms the simple Average forecast. This indicates that Apple’s log returns have no clear trend or seasonality, and the data behave almost like white noise, making the historical average the most reliable simple predictor.

d)

Why does ETS not work well for log return data? - Exponential smoothing models are unsuitable for log return data because returns have no clear level, trend, or seasonality — they behave almost like white noise. These models assume past patterns can predict future values, but financial returns are largely random, so ETS models offer no real forecasting advantage here.

h <- length(apple_test)

ets_specs <- list(
  "ETS(A,N,N)"  = list(model = "ANN", damped = FALSE), # Simple Exponential Smoothing
  "ETS(A,A,N)"  = list(model = "AAN", damped = FALSE), # Holt (trend)
  "ETS(A,Ad,N)" = list(model = "AAN", damped = TRUE)   # Holt-damped trend
)

fit_one <- function(lbl, spec) {
  fit <- ets(apple_train, model = spec$model, damped = spec$damped,
             allow.multiplicative.trend = TRUE)
  fc <- forecast(fit, h = h)
  a  <- accuracy(fc, apple_test)
  tibble(
    method = lbl,
    RMSE   = as.numeric(a["Test set","RMSE"]),
    MAE    = as.numeric(a["Test set","MAE"])
  )
}

ets_results <- bind_rows(
  lapply(names(ets_specs), function(lbl) fit_one(lbl, ets_specs[[lbl]]))
)


# Legg til ETS-radene
results <- bind_rows(results, ets_results)
print(results)
##           method       RMSE        MAE
## 1        Average 0.07162500 0.05901900
## 2          Drift 0.11407800 0.09504700
## 3          Naive 0.11406100 0.09503000
## 4 Seasonal Naive 0.08996300 0.07270400
## 5     ETS(A,N,N) 0.07162642 0.05901962
## 6     ETS(A,A,N) 0.07588771 0.06289294
## 7    ETS(A,Ad,N) 0.07579829 0.06278565

method RMSE MAE

1 Average 0.07161900 0.05899800 2 Drift 0.11405900 0.09502600 3 Naive 0.11404200 0.09500900 4 Seasonal Naive 0.08994000 0.07268300 5 ETS(A,N,N) 0.07161999 0.05899879 6 ETS(A,A,N) 0.07587665 0.06287215 7 ETS(A,Ad,N) 0.07578747 0.06276480

Comments: The exponential smoothing models (ETS) perform very similarly to the baseline methods, with ETS(A,N,N) producing nearly identical RMSE and MAE as the Average method. This shows that exponential smoothing provides no meaningful improvement for log return data. Because log returns lack a clear level, trend, or seasonality, the smoothing process has little structure to model. Models like ETS(A,A,N) and ETS(A,Ad,N) slightly worsen performance, as they introduce unnecessary trend components that do not exist in the data.

# Plot best ETS model
best_lbl <- ets_results |> arrange(RMSE) |> slice(1) |> pull(method)
best_spec <- ets_specs[[best_lbl]]
best_fit <- ets(apple_train, model = best_spec$model, damped = best_spec$damped,
                allow.multiplicative.trend = TRUE)
best_fc  <- forecast(best_fit, h = h)

autoplot(best_fc) +
  autolayer(apple_test, series = "Test") +
  ggtitle(paste0("Best ETS by RMSE: ", best_lbl)) +
  xlab("År") + ylab("Log-return")

e)

Based on the plot in Exercise 2, the Arima k,d,p:

suppressPackageStartupMessages({library(forecast); library(tibble); library(dplyr); library(ggplot2)})

h <- length(apple_test)

# Chosen from ACF/PACF of log returns: stationary, weak short memory -> AR(1)
p_manual <- 4  # AR
d_manual <- 0  # no differencing for returns
q_manual <- 0  # no MA term initially

fit_manual <- Arima(apple_train, order = c(p_manual, d_manual, q_manual), include.mean = TRUE)
fc_manual  <- forecast(fit_manual, h = h)

# Score and append to results
get_rmse_mae <- function(fc, test){
  a <- accuracy(fc, test)
  c(RMSE = as.numeric(a["Test set","RMSE"]),
    MAE  = as.numeric(a["Test set","MAE"]))
}
manual_row <- tibble(
  method = paste0("ARIMA(", p_manual, ",", d_manual, ",", q_manual, ")"),
  RMSE = get_rmse_mae(fc_manual, apple_test)["RMSE"],
  MAE  = get_rmse_mae(fc_manual, apple_test)["MAE"]
)
if (!exists("results")) results <- tibble(method=character(), RMSE=double(), MAE=double())
results <- bind_rows(results, manual_row)
print(results)
##           method       RMSE        MAE
## 1        Average 0.07162500 0.05901900
## 2          Drift 0.11407800 0.09504700
## 3          Naive 0.11406100 0.09503000
## 4 Seasonal Naive 0.08996300 0.07270400
## 5     ETS(A,N,N) 0.07162642 0.05901962
## 6     ETS(A,A,N) 0.07588771 0.06289294
## 7    ETS(A,Ad,N) 0.07579829 0.06278565
## 8   ARIMA(4,0,0) 0.07155605 0.05908402

Why 4,0,0?
- From the ACF and PACF plots of Apple’s log returns, the ACF shows no strong or slowly decaying correlations, indicating the series is stationary, so we set d = 0. - The PACF has a few small but noticeable spikes up to around lag 4, while the ACF cuts off quickly after lag 1. This pattern is typical of an autoregressive process of order 4, meaning the current value depends on the last four observations. - There is no clear moving-average structure (no tailing ACF), so we set q = 0

The autor arima model:

fit_auto <- auto.arima(apple_train, seasonal = FALSE, stepwise = TRUE, approximation = FALSE)
fc_auto  <- forecast(fit_auto, h = h)

auto_row <- tibble(
  method = paste0("auto.arima: ARIMA(", arimaorder(fit_auto)[1], ",", arimaorder(fit_auto)[2], ",", arimaorder(fit_auto)[3], ")"),
  RMSE = accuracy(fc_auto, apple_test)["Test set","RMSE"] |> as.numeric(),
  MAE  = accuracy(fc_auto, apple_test)["Test set","MAE"]  |> as.numeric()
)
results <- bind_rows(results, auto_row)
print(results)
##                     method       RMSE        MAE
## 1                  Average 0.07162500 0.05901900
## 2                    Drift 0.11407800 0.09504700
## 3                    Naive 0.11406100 0.09503000
## 4           Seasonal Naive 0.08996300 0.07270400
## 5               ETS(A,N,N) 0.07162642 0.05901962
## 6               ETS(A,A,N) 0.07588771 0.06289294
## 7              ETS(A,Ad,N) 0.07579829 0.06278565
## 8             ARIMA(4,0,0) 0.07155605 0.05908402
## 9 auto.arima: ARIMA(0,0,0) 0.07162501 0.05901876
cat("Manual order: ARIMA(", p_manual, ",", d_manual, ",", q_manual, ")\n", sep = "")
## Manual order: ARIMA(4,0,0)
cat("Auto order:   ARIMA(", paste(arimaorder(fit_auto)[c("p","d","q")], collapse=","), ")\n", sep = "")
## Auto order:   ARIMA(0,0,0)

Plot:

p1 <- autoplot(fc_manual) + autolayer(apple_test, series = "Test") +
  ggtitle(paste0("Manual ARIMA(", p_manual, ",", d_manual, ",", q_manual, ")")) +
  xlab("Year") + ylab("Log return")

print(p1)

Comments on Arima(4,0,0) plot: - The forecast line is nearly flat around zero meaning that the model predicts that the future return will stay close to the long term average. - The test-data lie inside the confidence band, meaning the model’s uncertainty estimate is reasonable. - The forecast itself does not track the short-term movements in the test data.

p2 <- autoplot(fc_auto) + autolayer(apple_test, series = "Test") +
  ggtitle(paste0("Auto ARIMA: ARIMA(", paste(arimaorder(fit_auto)[c("p","d","q")], collapse=","), ")")) +
  xlab("Year") + ylab("Log return")

print(p2)

Comments on Auto Arima plot: - The forecast line is also nearly flat around zero, indicating that the model predicts future returns will remain close to the long-term average. - The test data points mostly lie within the confidence intervals, suggesting that the model’s uncertainty estimates are appropriate. - Similar to the manual ARIMA model, the forecast does not capture the short-term fluctuations in the test data.

RMSE and MAE on test data:

cat("\n-- Residual diagnostics: manual --\n")
## 
## -- Residual diagnostics: manual --
print(Box.test(residuals(fit_manual), lag = 12, type = "Ljung-Box", fitdf = p_manual + q_manual))
## 
##  Box-Ljung test
## 
## data:  residuals(fit_manual)
## X-squared = 15.383, df = 8, p-value = 0.05211
checkresiduals(fit_manual)

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(4,0,0) with non-zero mean
## Q* = 23.014, df = 20, p-value = 0.2881
## 
## Model df: 4.   Total lags used: 24

Box-Ljung test data: residuals(fit_manual) X-squared = 14.981, df = 8, p-value = 0.05952

Ljung-Box test Q* = 22.809, df = 20, p-value = 0.2983 Model df: 4. Total lags used: 24

Comments: - The p-values (0.059 and 0.298) are both above 0.05, meaning we fail to reject the null hypothesis of no autocorrelation. - The residuals appear to be white noise, indicating that the ARIMA(4,0,0) model adequately captures the structure in the data. - No significant patterns remain unexplained by the model.

cat("\n-- Residual diagnostics: auto --\n")
## 
## -- Residual diagnostics: auto --
fitdf_auto <- arimaorder(fit_auto)["p"] + arimaorder(fit_auto)["q"]
print(Box.test(residuals(fit_auto), lag = 12, type = "Ljung-Box", fitdf = fitdf_auto))
## 
##  Box-Ljung test
## 
## data:  residuals(fit_auto)
## X-squared = 19.826, df.p = 12, p-value = 0.07045
checkresiduals(fit_auto)

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(0,0,0) with non-zero mean
## Q* = 28, df = 24, p-value = 0.26
## 
## Model df: 0.   Total lags used: 24

Box-Ljung test data: residuals(fit_auto) X-squared = 19.655, df.p = 12, p-value = 0.0739 >

Ljung-Box test data: Residuals from ARIMA(0,0,0) Q* = 28.018, df = 24, p-value = 0.2593

Comments: - Both p-values (0.074 and 0.259) are well above 0.05, so we fail to reject the null hypothesis of no autocorrelation. - The residuals behave like white noise, meaning the ARIMA(0,0,0) model adequately fits the stationary log-return data. - No significant autocorrelation remains, confirming that the returns are essentially random around their mean.

RMSE and MAE results: For both models: Top plot-time series of residuals: - The residuals appear to fluctuate randomly around zero, indicating no obvious patterns or trends, meaning the model doesn’t systematically over- or under-predict.

ACF plot: - The ACF of the residuals shows that all autocorrelation values are within the confidence bounds, suggesting no significant autocorrelation remains. This indicates that the model has effectively captured the temporal dependencies in the data.

Bottom-right histogram: - The histogram of the residuals appears approximately normally distributed, which is a common assumption in time series modeling. This suggests that the model’s errors are random and not biased.

f)

results_sorted <- results %>%
  arrange(RMSE)

print(results_sorted)
##                     method       RMSE        MAE
## 1             ARIMA(4,0,0) 0.07155605 0.05908402
## 2                  Average 0.07162500 0.05901900
## 3 auto.arima: ARIMA(0,0,0) 0.07162501 0.05901876
## 4               ETS(A,N,N) 0.07162642 0.05901962
## 5              ETS(A,Ad,N) 0.07579829 0.06278565
## 6               ETS(A,A,N) 0.07588771 0.06289294
## 7           Seasonal Naive 0.08996300 0.07270400
## 8                    Naive 0.11406100 0.09503000
## 9                    Drift 0.11407800 0.09504700

Results:
method RMSE MAE ———————————————— 1 ARIMA(4,0,0) 0.07154955 0.05906318 2 auto.arima: ARIMA(0,0,0) 0.07161857 0.05899792 3 Average 0.07161900 0.05899800 4 ETS(A,N,N) 0.07161999 0.05899879 5 ETS(A,Ad,N) 0.07578747 0.06276480 6 ETS(A,A,N) 0.07587665 0.06287215 7 Seasonal Naive 0.08994000 0.07268300 8 Naive 0.11404200 0.09500900 9 Drift 0.11405900 0.09502600

Comments: It looks like the ETS model edges out the average method, but the improvement is negligible, all models are effectively forecasting random noise around zero, so it’s inherently difficult to beat a simple baseline.

Exercise 4

a)

symbols <- c("^GSPC", "^VIX")

sp500_vix_daily <- get_finance_data(symbols, "2000-01-01", Sys.Date())

# Check structure
glimpse(sp500_vix_daily)
## Rows: 13,004
## Columns: 5
## $ date     <date> 2000-01-03, 2000-01-04, 2000-01-05, 2000-01-06, 2000-01-07, …
## $ symbol   <fct> ^GSPC, ^GSPC, ^GSPC, ^GSPC, ^GSPC, ^GSPC, ^GSPC, ^GSPC, ^GSPC…
## $ open     <dbl> 1469.25, 1455.22, 1399.42, 1402.11, 1403.45, 1441.47, 1457.60…
## $ low      <dbl> 1438.36, 1397.43, 1377.68, 1392.10, 1400.73, 1441.47, 1434.42…
## $ adjusted <dbl> 1455.22, 1399.42, 1402.11, 1403.45, 1441.47, 1457.60, 1438.56…

Convert daily data to monthly data by taking last trading day of each month

sp500_vix_monthly <- sp500_vix_daily %>%
  dplyr::mutate(year_month = lubridate::floor_date(date, "month")) %>%  
  dplyr::group_by(symbol, year_month) %>%
  dplyr::slice_max(order_by = date, n = 1, with_ties = FALSE) %>%      
  dplyr::ungroup() %>%
  dplyr::select(date, symbol, adjusted) %>%
  dplyr::arrange(symbol, date)

head(sp500_vix_monthly)
## # A tibble: 6 × 3
##   date       symbol adjusted
##   <date>     <fct>     <dbl>
## 1 2000-01-31 ^GSPC     1394.
## 2 2000-02-29 ^GSPC     1366.
## 3 2000-03-31 ^GSPC     1499.
## 4 2000-04-28 ^GSPC     1452.
## 5 2000-05-31 ^GSPC     1421.
## 6 2000-06-30 ^GSPC     1455.
# Pivot to wide format for easier analysis
sp500_vix_wide <- sp500_vix_monthly %>%
  pivot_wider(names_from = symbol, values_from = adjusted) %>%
  arrange(date)

head(sp500_vix_wide)
## # A tibble: 6 × 3
##   date       `^GSPC` `^VIX`
##   <date>       <dbl>  <dbl>
## 1 2000-01-31   1394.   25.0
## 2 2000-02-29   1366.   23.4
## 3 2000-03-31   1499.   24.1
## 4 2000-04-28   1452.   26.2
## 5 2000-05-31   1421.   23.6
## 6 2000-06-30   1455.   19.5
# Check how many missing values in each column
colSums(is.na(sp500_vix_wide))
##  date ^GSPC  ^VIX 
##     0     1     1
# Fill missing values using last observation carried forward (LOCF) due to non-trading days
sp500_vix_wide <- sp500_vix_wide %>%
  mutate(across(everything(), ~na.locf(.x, na.rm = FALSE)))

# Confirm that NAs are gone
colSums(is.na(sp500_vix_wide))
##  date ^GSPC  ^VIX 
##     0     0     0

Summary statistics:

summary(sp500_vix_wide[, -1])
##      ^GSPC             ^VIX      
##  Min.   : 735.1   Min.   : 9.51  
##  1st Qu.:1217.1   1st Qu.:14.04  
##  Median :1528.7   Median :17.79  
##  Mean   :2275.7   Mean   :19.88  
##  3rd Qu.:2917.1   3rd Qu.:23.63  
##  Max.   :6840.2   Max.   :59.89
 ^GSPC             ^VIX

Min. : 735.1 Min. : 9.51
1st Qu.:1217.1 1st Qu.:14.04
Median :1528.7 Median :17.79
Mean :2276.6 Mean :19.88
3rd Qu.:2917.1 3rd Qu.:23.63
Max. :6852.0 Max. :59.89

Comments: - The S&P 500 (ˆGSPC) index ranges from about 735 to 6852, showing strong long-term growth since 2000. - The VIX (volatility index) ranges from 9.5 to nearly 60, indicating large swings in market uncertainty. - The median VIX (≈18) suggests that typical market volatility is moderate, while the high maximum values reflect short periods of financial stress (e.g., 2008 crisis or 2020 pandemic). - Overall, the S&P 500 trend is upward, while the VIX captures episodes of fear and market turbulence.

ggplot(sp500_vix_wide, aes(x = date)) +
  geom_line(aes(y = `^GSPC`, color = "S&P500")) +
  geom_line(aes(y = `^VIX` * 100, color = "VIX × 100 (scaled)"), linetype = "dashed") +
  scale_y_continuous(
    name = "S&P 500 (index level)",
    sec.axis = sec_axis(~./100, name = "VIX (volatility index)")
  ) +
  labs(
    title = "S&P 500 and VIX Monthly Closing Prices (2000–Today)",
    x = "Year"
  ) +
  theme_minimal() +
  theme(
    plot.title = element_text(size = 16, hjust = 0.5),
    axis.title = element_text(size = 13),
    legend.title = element_blank(),
    legend.position = "bottom"
  )

Comments on the plot: - The S&P 500 shows a clear upward trend over time, reflecting long-term economic growth. - The VIX, which measures market volatility and investor fear, moves inversely to the S&P 500. - Periods of market downturns (e.g., 2008–2009 financial crisis and 2020 COVID-19 crash) correspond to sharp spikes in the VIX, indicating increased uncertainty. - In stable bull markets, the VIX remains low, while the S&P 500 steadily rises. - Overall, the plot highlights the negative relationship between stock market performance and volatility: when fear rises, prices fall.

b)

Compute log returns/differences:

names(sp500_vix_wide)
## [1] "date"  "^GSPC" "^VIX"
sp500_vix_log <- sp500_vix_wide %>%
  dplyr::arrange(date) %>%
  dplyr::mutate(
    log_return_sp500 = log(`^GSPC` / dplyr::lag(`^GSPC`)),
    log_diff_vix     = log(`^VIX`  / dplyr::lag(`^VIX`))
  ) %>%
  dplyr::select(date, log_return_sp500, log_diff_vix) %>%
  tidyr::drop_na()

head(sp500_vix_log)
## # A tibble: 6 × 3
##   date       log_return_sp500 log_diff_vix
##   <date>                <dbl>        <dbl>
## 1 2000-02-29          -0.0203      -0.0654
## 2 2000-03-31           0.0923       0.0312
## 3 2000-04-28          -0.0313       0.0831
## 4 2000-05-31          -0.0222      -0.102 
## 5 2000-06-30           0.0237      -0.191 
## 6 2000-07-31          -0.0165       0.0596
apple_log <- apple_log_returns %>%
  dplyr::rename(log_return_apple = log_return) %>%
  dplyr::select(date, log_return_apple)

head(apple_log)
## # A tibble: 6 × 2
##   date       log_return_apple
##   <date>                <dbl>
## 1 2000-02-29           0.0997
## 2 2000-03-31           0.170 
## 3 2000-04-28          -0.0905
## 4 2000-05-31          -0.390 
## 5 2000-06-30           0.221 
## 6 2000-07-31          -0.0303

Merge all three series

combined_returns <- sp500_vix_log %>%
  left_join(apple_log, by = "date")

head(combined_returns)
## # A tibble: 6 × 4
##   date       log_return_sp500 log_diff_vix log_return_apple
##   <date>                <dbl>        <dbl>            <dbl>
## 1 2000-02-29          -0.0203      -0.0654           0.0997
## 2 2000-03-31           0.0923       0.0312           0.170 
## 3 2000-04-28          -0.0313       0.0831          -0.0905
## 4 2000-05-31          -0.0222      -0.102           -0.390 
## 5 2000-06-30           0.0237      -0.191            0.221 
## 6 2000-07-31          -0.0165       0.0596          -0.0303

Visualise log returns

library(scales)

ggplot(combined_returns, aes(x = date)) +
  geom_line(aes(y = log_return_sp500, color = "S&P 500")) +
  geom_line(aes(y = log_diff_vix, color = "VIX (log diff)"), linetype = "dashed") +
  geom_line(aes(y = log_return_apple, color = "Apple"), alpha = 0.7) +
  labs(
    title = "Monthly Log Returns – Apple, S&P 500, and VIX (2000–Today)",
    x = "Year",
    y = "Log return / difference"
  ) +
  scale_color_manual(values = c("Apple" = "#1f77b4", "S&P 500" = "#e41a1c", "VIX (log diff)" = "#4daf4a")) +
  scale_y_continuous(labels = percent_format(accuracy = 1)) +
  theme_minimal() +
  theme(
    plot.title = element_text(size = 15, hjust = 0.5),
    legend.title = element_blank(),
    legend.position = "bottom"
  )

Zoom in on periode:

# Example: zoom in on 2007–2010 to show the financial crisis
ggplot(filter(combined_returns, date >= as.Date("2007-01-01"), date <= as.Date("2010-12-31")),
       aes(x = date)) +
  geom_line(aes(y = log_return_sp500, color = "S&P 500")) +
  geom_line(aes(y = log_diff_vix, color = "VIX (log diff)"), linetype = "dashed") +
  geom_line(aes(y = log_return_apple, color = "Apple"), alpha = 0.7) +
  labs(title = "Zoom: Log Returns during 2008 Financial Crisis",
       x = "Year", y = "Log return / difference") +
  scale_color_manual(values = c("Apple" = "#1f77b4", "S&P 500" = "#e41a1c", "VIX (log diff)" = "#4daf4a")) +
  theme_minimal() +
  theme(plot.title = element_text(size = 14, hjust = 0.5),
        legend.title = element_blank(),
        legend.position = "bottom")

Comments on plot: - The S&P 500 and Apple show very similar movements, confirming that Apple’s stock follows overall market trends. - Apple’s returns are more volatile, with larger positive and negative swings, reflecting higher risk. - The VIX tends to spike exactly when both Apple and the S&P 500 show large negative returns, especially during crisis periods (2008–2009 and early 2020). - This demonstrates the classic negative correlation between market returns and volatility: when fear and uncertainty rise, stock prices fall. - Outside of crisis periods, all three series fluctuate around zero, suggesting returns are largely random with short-lived volatility shocks.

c)

Prepare each series with year/month keys:

# S&P500 + VIX monthly (from 4a)
sp500_vix_log <- sp500_vix_wide %>%
  dplyr::arrange(date) %>%
  dplyr::mutate(
    log_return_sp500 = log(`^GSPC` / dplyr::lag(`^GSPC`)),
    log_diff_vix     = log(`^VIX`  / dplyr::lag(`^VIX`)),
    year  = lubridate::year(date),
    month = lubridate::month(date)
  ) %>%
  dplyr::select(date, year, month, log_return_sp500, log_diff_vix) %>%
  tidyr::drop_na()

# Company (Apple) monthly log returns (pick your available source)
if (exists("apple_log_returns")) {
  company_log <- apple_log_returns %>%
    dplyr::rename(log_return_company = log_return) %>%
    dplyr::mutate(year = lubridate::year(date),
                  month = lubridate::month(date)) %>%
    dplyr::select(date, year, month, log_return_company)
} else if (exists("log_returns_wide") && "AAPL" %in% names(log_returns_wide)) {
  company_log <- log_returns_wide %>%
    dplyr::rename(log_return_company = AAPL) %>%
    dplyr::mutate(year = lubridate::year(date),
                  month = lubridate::month(date)) %>%
    dplyr::select(date, year, month, log_return_company)
} else if (exists("apple_returns_ts")) {
  company_log <- tibble::tibble(
    date  = as.Date(zoo::as.yearmon(time(apple_returns_ts))),
    log_return_company = as.numeric(apple_returns_ts)
  ) %>%
    dplyr::mutate(year = lubridate::year(date),
                  month = lubridate::month(date)) %>%
    dplyr::select(date, year, month, log_return_company)
} else {
  stop("Could not find the company monthly log returns.")
}

Inner join on year/month then build a multivariate ts:

mv_monthly <- sp500_vix_log %>%
  dplyr::inner_join(company_log, by = c("year","month")) %>%
  dplyr::arrange(date.x) %>%
  dplyr::transmute(
    date = date.x,
    log_return_company,
    log_return_sp500,
    log_diff_vix
  )

# Create an mts (multivariate ts) with frequency = 12
start_yr  <- lubridate::year(min(mv_monthly$date))
start_mon <- lubridate::month(min(mv_monthly$date))
X_mts <- ts(as.matrix(mv_monthly[, c("log_return_company","log_return_sp500","log_diff_vix")]),
            start = c(start_yr, start_mon), frequency = 12)
colnames(X_mts) <- c("company","sp500","vix_logdiff")

Split into train/test (last 48 months = test)

h <- 48
n <- nrow(X_mts)
stopifnot(n > h)

ts_time <- time(X_mts)
X_train <- window(X_mts, end   = ts_time[n - h])
X_test  <- window(X_mts, start = ts_time[n - h + 1])

# Quick sanity checks
dim(X_train); dim(X_test)   # X_test should have 48 rows
## [1] 263   3
## [1] 48  3
head(X_train); head(X_test)
##              company       sp500 vix_logdiff
## Feb 2000  0.09968247 -0.02031300 -0.06542067
## Mar 2000  0.16960942  0.09232375  0.03117353
## Apr 2000 -0.09049047 -0.03127991  0.08313272
## May 2000 -0.38996851 -0.02215875 -0.10239634
## Jun 2000  0.22075967  0.02365163 -0.19089941
## Jul 2000 -0.03028710 -0.01647627  0.05960050
##              company         sp500 vix_logdiff
## Jan 2022 -0.01583675 -5.401823e-02  0.36598114
## Feb 2022 -0.05558245 -3.186276e-02  0.19413231
## Mar 2022  0.05588259  3.514829e-02 -0.38283750
## Apr 2022 -0.10217763 -9.206782e-02  0.48520853
## May 2022 -0.05603722  5.317629e-05 -0.24317827
## Jun 2022 -0.08493695 -8.765158e-02  0.09186778

d)

Fit the three ARIMAX models:

library(forecast)
library(dplyr)

# Response = company log return (Apple)
y_train <- X_train[, "company"]
y_test  <- X_test[,  "company"]

# Exogenous regressors
xreg_train_sp500 <- X_train[, "sp500", drop = FALSE]
xreg_test_sp500  <- X_test[,  "sp500",  drop = FALSE]

xreg_train_vix <- X_train[, "vix_logdiff", drop = FALSE]
xreg_test_vix  <- X_test[,  "vix_logdiff",  drop = FALSE]

xreg_train_both <- X_train[, c("sp500","vix_logdiff")]
xreg_test_both  <- X_test[,  c("sp500","vix_logdiff")]

# 1) ARIMAX with S&P500
fit_sp500 <- auto.arima(y_train, xreg = xreg_train_sp500, seasonal = FALSE)
fc_sp500  <- forecast(fit_sp500, xreg = xreg_test_sp500, h = nrow(xreg_test_sp500))

# 2) ARIMAX with VIX
fit_vix <- auto.arima(y_train, xreg = xreg_train_vix, seasonal = FALSE)
fc_vix  <- forecast(fit_vix, xreg = xreg_test_vix, h = nrow(xreg_test_vix))

# 3) ARIMAX with both S&P500 + VIX
fit_both <- auto.arima(y_train, xreg = xreg_train_both, seasonal = FALSE)
fc_both  <- forecast(fit_both, xreg = xreg_test_both, h = nrow(xreg_test_both))

Compute RMSE & MAE and append to results

get_rmse_mae <- function(fc, actual) {
  data.frame(
    RMSE = sqrt(mean((actual - fc$mean)^2, na.rm = TRUE)),
    MAE  = mean(abs(actual - fc$mean), na.rm = TRUE)
  )
}

arimax_results <- bind_rows(
  cbind(method = "ARIMAX (S&P500)",       get_rmse_mae(fc_sp500, y_test)),
  cbind(method = "ARIMAX (VIX)",          get_rmse_mae(fc_vix,   y_test)),
  cbind(method = "ARIMAX (S&P500 + VIX)", get_rmse_mae(fc_both,  y_test))
)

# Append to the main results from Exercise 3
if (!exists("results")) results <- arimax_results else results <- bind_rows(results, arimax_results)

results <- distinct(results)
print(results %>% arrange(RMSE))
##                      method       RMSE        MAE
## 1     ARIMAX (S&P500 + VIX) 0.05319619 0.04145037
## 2           ARIMAX (S&P500) 0.05339243 0.04191533
## 3              ARIMAX (VIX) 0.06345857 0.05239851
## 4              ARIMA(4,0,0) 0.07155605 0.05908402
## 5                   Average 0.07162500 0.05901900
## 6  auto.arima: ARIMA(0,0,0) 0.07162501 0.05901876
## 7                ETS(A,N,N) 0.07162642 0.05901962
## 8               ETS(A,Ad,N) 0.07579829 0.06278565
## 9                ETS(A,A,N) 0.07588771 0.06289294
## 10           Seasonal Naive 0.08996300 0.07270400
## 11                    Naive 0.11406100 0.09503000
## 12                    Drift 0.11407800 0.09504700

Visualize forecasts (vs actual returns) ARIMAX with S&P500

autoplot(fc_sp500) + autolayer(y_test, series="Actual") +
  ggtitle("ARIMAX Forecast – Exogenous: S&P500") +
  xlab("Time") + ylab("Company log return")

autoplot(fc_vix) + autolayer(y_test, series="Actual") +
  ggtitle("ARIMAX Forecast – Exogenous: VIX") +
  xlab("Time") + ylab("Company log return")

autoplot(fc_both) + autolayer(y_test, series="Actual") +
  ggtitle("ARIMAX Forecast – Exogenous: S&P500 + VIX") +
  xlab("Time") + ylab("Company log return")

Comments on the plots: - The ARIMAX model with S&P 500 and VIX clearly outperforms the baseline forecasts. - It captures short-term movements in Apple’s returns more accurately, as the predicted values follow the actual data closely. - Including market and volatility information improves the model’s explanatory power, even though forecast uncertainty remains high due to the inherent randomness of financial returns.

Residual diagnostics

checkresiduals(fit_sp500)

## 
##  Ljung-Box test
## 
## data:  Residuals from Regression with ARIMA(0,0,0) errors
## Q* = 23.762, df = 24, p-value = 0.4753
## 
## Model df: 0.   Total lags used: 24
checkresiduals(fit_vix)

## 
##  Ljung-Box test
## 
## data:  Residuals from Regression with ARIMA(0,0,0) errors
## Q* = 25.875, df = 24, p-value = 0.3595
## 
## Model df: 0.   Total lags used: 24
checkresiduals(fit_both)

## 
##  Ljung-Box test
## 
## data:  Residuals from Regression with ARIMA(0,0,0) errors
## Q* = 23.773, df = 24, p-value = 0.4746
## 
## Model df: 0.   Total lags used: 24

Comments on the residual: Top plot (residuals over time): - The residuals fluctuate randomly around zero with no visible trend or clustering. This indicates that the model has captured the main dynamics in the data and that no systematic patterns remain.

Bottom-left (ACF of residuals): - Most autocorrelation values lie within the blue confidence bounds, suggesting no significant autocorrelation. The residuals behave like white noise, fulfilling the ARIMAX model assumptions.

Bottom-right (histogram + density): - The residuals are approximately normally distributed, centered near zero. The shape is slightly leptokurtic (heavy tails), which is typical for financial return data.

e)

# Combine everything (baseline + ETS + ARIMA + ARIMAX)
results_sorted <- results %>%
  dplyr::arrange(RMSE)

print(results_sorted)
##                      method       RMSE        MAE
## 1     ARIMAX (S&P500 + VIX) 0.05319619 0.04145037
## 2           ARIMAX (S&P500) 0.05339243 0.04191533
## 3              ARIMAX (VIX) 0.06345857 0.05239851
## 4              ARIMA(4,0,0) 0.07155605 0.05908402
## 5                   Average 0.07162500 0.05901900
## 6  auto.arima: ARIMA(0,0,0) 0.07162501 0.05901876
## 7                ETS(A,N,N) 0.07162642 0.05901962
## 8               ETS(A,Ad,N) 0.07579829 0.06278565
## 9                ETS(A,A,N) 0.07588771 0.06289294
## 10           Seasonal Naive 0.08996300 0.07270400
## 11                    Naive 0.11406100 0.09503000
## 12                    Drift 0.11407800 0.09504700

Comments: The ARIMAX models clearly outperform all baseline, ETS, and ARIMA models in terms of both RMSE and MAE. The ARIMAX model with S&P 500 and VIX provides the best overall performance, showing that incorporating market and volatility information improves forecasts of Apple’s returns. The S&P 500 alone also adds strong predictive power, reflecting Apple’s close correlation with the broader market. In contrast, models relying only on past Apple data (ARIMA/ETS/baseline) perform worse and mostly predict around the mean.

Pros of exogenous variables: - Capture relationships between Apple’s returns, the overall market, and volatility. - Improve short-term forecast accuracy and responsiveness to market changes. - Provide more realistic forecasts during periods of high volatility.

Cons: - The improvement, while significant, remains moderate — returns are still largely random. - More complex models can overfit or rely on variables that are themselves uncertain.

Summary: Yes — the ARIMAX models, especially the one including both S&P 500 and VIX, successfully beat the baseline models from Exercise 3. Adding exogenous variables yields meaningful performance gains, but forecasting financial returns remains inherently difficult due to their noisy and unpredictable nature.